#ifndef _contest_hungarian_h_
#define _contest_hungarian_h_

#include<contest/vectord.h>

namespace contest
{

#define INF 2147483647

/**
 * The hungarian algorithm.
 * Find a permutation which sums is minimum.
 **/
class hungarian
{
private:
    size_t n;
    vectord<int32_t, 2> aa;
    vector< int32_t > u, v;
    vector< size_t > col_mate, row_mate, parent;

    pair<size_t, size_t> augment()
    {
        size_t FREE = n + 1;
        vector< size_t > unchosen(n);
        vector< size_t > slack_row(n);
        vector< int32_t > slack(n, INF);
        for (size_t c=0; c<n; ++c) 
            parent[c] = FREE;

        size_t t = 0;
        for (size_t r=0; r<n; ++r) 
            if (col_mate[r] == FREE) 
                unchosen[t++] = r;

        while (true) 
        {
            for (size_t q=0; q<t; ++q)
            {
                size_t r = unchosen[q];
                for (size_t c=0; c<n; ++c)
                    if (slack[c] != 0 && aa(r, c) - u[r] + v[c] < slack[c]) 
                    {
                        slack[c] = aa(r, c) - u[r] + v[c];
                        if (slack[c] == 0) 
                        {
                            if (row_mate[c] == FREE)
                            {
                                return pair<size_t, size_t>(r, c);
                            }
                            parent[c] = r;
                            unchosen[t++] = row_mate[c];
                        } 
                        else slack_row[c] = r;
                    }
            }
            int32_t s = INF;
            for (size_t c=0; c<n; ++c)
                if (slack[c] != 0 && slack[c] < s) 
                    s = slack[c];
            for (size_t q=0; q<t; ++q) 
                u[unchosen[q]]+=s;
            for (size_t c=0; c<n; ++c)
                if (slack[c] != 0) 
                {
                    slack[c] -= s;
                    if (slack[c] == 0) 
                    {
                        size_t r = slack_row[c];
                        if (row_mate[c] == FREE) 
                        {
                            for (size_t j=c+1; j<n; ++j )
                                if (slack[j] == 0) 
                                    v[j] += s;
                            return pair<size_t, size_t>(r, c);
                        }
                        else
                        {
                            parent[c] = r;
                            unchosen[t++] = row_mate[c];
                        }
                    }
                }//slack[l]!=0
                else v[c] += s; //increasing zeros
        }//while (1)
    }

public:
    hungarian(vectord<int32_t, 2>& mtx, const size_t side):
        aa(mtx), n(side), 
        u(n, 0), v(n, 0),
        col_mate(n, n + 1), row_mate(n, n + 1), parent(n, n + 1)
    {
        size_t FREE = n + 1;

        size_t unmatched = 0;
        for (size_t r=0; r<n; ++r) 
        {
            size_t pos = 0;
            for (size_t c=1; c<n; ++c)
                if (aa(r, c) < aa(r, pos)) 
                    pos = c;
            u[r] = aa(r, pos);
            if (row_mate[pos] == FREE) 
            {
                row_mate[pos] = r;
                col_mate[r] = pos;
            } 
            else ++unmatched;
        }

        if (unmatched == 0) return;

        while (true) 
        {
            pair<size_t, size_t> tuple = augment();
            size_t row = tuple.first;
            size_t col = tuple.second;
            while (true) 
            {
                size_t j = col_mate[row];
                col_mate[row] = col;
                row_mate[col] = row;
                if (j == FREE) break;
                row = parent[j];
                col = j;
            }
            if (--unmatched == 0) break;
        }
    }

    const size_t get_column(const size_t row) const
    {
        return col_mate[row];
    }
};

}

#endif
